Diesel Cycle in Mathematica for Education

David ((Smrčka)^1), Zdeněk (Pucholt^( 2))

Palacký University in Olomouc, Czech Republic

http://www.prf.upol.cz/en/

david.smrcka@upol.(Subscript[cz, 1]), zdenek.pucholt@upol.(Subscript[cz, 2])

This presentation is supported by project
Moderní přístup k aplikaci matematických dovedností v přírodovědných a ekonomických oborech MATHSCIENCE
CZ.1.07/2.2.00/28.0168

(1) Introduction to Diesel cycle

thermodynamic cycle is approximating the pressure and volume of the combustion chamber of the diesel engine, invented by Rudolph Diesel in 1893

Rudolph Diesel

pV diagram consists of 4 phases:

process 1 to 2 is isentropic compression of the fluid

(Subscript[W, in]) done by piston compressing working fluid

process 2 to 3 is reversible constant pressure heating

(Subscript[Q, in]) by the combustion of the fuel

process 3 to 4 is isentropic expansion

(Subscript[W, out]) by the working fluid expanding on to the piston → produces usable torque

process 4 to 1 is reversible constant volume cooling

(Subscript[Q, out]) by venting the air

pV diagram of Diesel cycle

The maximum (ideal) thermal efficiency (Subscript[η, max]) of Diesel cycle depends on:

cut-off ratio α=(Subscript[V, 3])/(Subscript[V, 2])

compression ration β=(Subscript[V, 1])/(Subscript[V, 2])

ratio of specific heats γ=(Subscript[C, p])/(Subscript[C, V])

(Subscript[η, max])=1-(1/(β^(γ-1)))((((α^γ)-1)/(γ(α-1))))

the actual thermal efficiency will be significantly lower due to heat and friction losses

(2) Diesel engine

Diesel engine is a heat engine: it converts heat into work

used in cars, power generation, diesel-electric locomotives, and submarines

highest thermal efficiency of any standard internal or external combustion engine due to its very high compression ratio and inherent lean burn which facilitates unburnt gasses to scavenge waste heat

low-speed diesel engines can have a thermal efficiency that exceeds 50 % (as used in ships and other applications where overall engine weight is relatively unimportant)

(3) Diesel engine ~ “briefly” history

1891: Herbert Akroyd Stuart invents the first internal combustion engine to use a pressurised fuel injection system

1892: February 23, Rudolf Diesel obtained a patent (RP 67207) titled “Arbeitsverfahren und Ausführungsart für Verbrennungsmaschinen”

10. 8. 1893: Date of invention Diesel engine - 121 years ago (source: WolframAlpha)

1897: Rudolf Diesel built his well-known high-compression prototype engine in 1897

Diesel prototype engine in 1897

1903: Two first diesel-powered ships were launched, both for river and canal operations: Petite-Pierre in France, powered by Dyckhoff-built diesels, and Vandal tanker in Russia

1904: The French built the first diesel submarine, the Z.

1905: Four diesel engine turbochargers and intercoolers were manufactured

1930: First US diesel-power passenger car (Cummins powered Packard) built in Columbus, Indiana (USA)

1933: First European passenger cars with diesel engines (Citroën Rosalie)

Manipulate[
GraphicsRow[
{dieselPV[p1,V1,β,α,γ],dieselPT[p1,V1,T1,β,α,γ],dieselTS[p1,V1,T1,0,β,α,γ,1.005]},
ImageSize->1000],

{{V1,875,"Initial volume [(cm^3)]"},10},{{p1,100,"Initial pressure [kPa]"},100},{{T1,293,"Initial temperature [K]"},293},{γ,1.4,ControlType->None},{{β,5,"Compression ratio"},.1,15,.1},{{α,1.5,"Cut-off ratio"},0,β-.1,.1},
Row[{Style["Efficiency of the cycle: ",20,Bold],Style[Dynamic[efficiency[α,β,γ]*100],20,Bold],Style[" %",20,Bold]}]
]

Manipulate[
Grid[{
{dieselPVTime[p1,V1,β,α,γ,t,350],dieselPTTime[p1,V1,T1,β,α,γ,t,350],dieselTSTime[p1,V1,T1,0,β,α,γ,1.005,t,350]},
{pistCiferniky[p1,V1,T1,β,α,γ,t,800],SpanFromLeft,SpanFromLeft}
}],

{{V1,875,"Initial volume [(cm^3)]"},10},{{p1,100,"Initial pressure [kPa]"},100},{{T1,293,"Initial temperature [K]"},293},{γ,1.4,ControlType->None},{{β,5,"Compression ratio"},.1,15,.1},{{α,1.5,"Cut-off ratio"},0,β-.1,.1},
{{t,0,"time"},0,40},
Row[{Style["Efficiency of the cycle: ",20,Bold],Style[Dynamic[efficiency[α,β,γ]*100],20,Bold],Style[" %",20,Bold]}]]

Manipulate[
Grid[{
{valecAhridel[p1,V1,β,α,γ,t],TabView[{"p-V diagram"->dieselPVTime[p1,V1,β,α,γ,t,500],"p-T diagram"->dieselPTTime[p1,V1,T1,β,α,γ,t,500],"T-s diagram"->dieselTSTime[p1,V1,T1,0,β,α,γ,1.005,t,500]}]},
{SpanFromAbove,pistCiferniky[p1,V1,T1,β,α,γ,t,600]}
}],

{{V1,875,"Initial volume [(cm^3)]"},10},{{p1,100,"Initial pressure [kPa]"},100},{{T1,293,"Initial temperature [K]"},293},{γ,1.4,ControlType->None},{{β,5,"Compression ratio"},.1,15,.1},{{α,1.5,"Cut-off ratio"},0,β-.1,.1},
{{t,0,"time"},0,40},
Row[{Style["Efficiency of the cycle: ",20,Bold],Style[Dynamic[efficiency[α,β,γ]*100],20,Bold],Style[" %",20,Bold]}]]

Thank you for your attention

Initialization

efficiency[α_,β_,γ_]:=1-(1/(β^(γ-1)))((((α^γ)-1)/(γ(α-1))))

dieselPV[p1_,V1_,ϵ_,ϕ_,κ_]:=
Module[{p2,p4,V2,V3,T2,T3,T4,pV1,pV2,pV3,pV4,adiKomp,adiExp,plot},
V2=V1/ϵ;V3=V2*ϕ;
p2=p1*(ϵ^κ);p4=p2*(V3/V1)^κ;

adiKomp=p1*(V1^κ);adiExp=p2*(V3^κ);
pV1={V1,p1};
pV2={V2,p2};
pV3={V3,p2};
pV4={V1,p4};
plot=Table[{V,adiKomp*V^(-κ)},{V,V2,V1,.001}];

Show[{
Plot[adiKomp*V^(-κ),{V,V1,V2},PlotStyle->{Blue,Thick}],
Plot[adiExp*V^(-κ),{V,V3,V1},PlotStyle->{Blue,Thick}],
Graphics[{
Blue,Thick,Line[{{pV2,pV3},{pV4,pV1}}]
}]
},AxesOrigin->{0,0},PlotRange->All,AspectRatio->1,PlotLabel->"p-V diagram"]
]

dieselPVTime[p1_,V1_,ϵ_,ϕ_,κ_,time_,imgS_]:=
Module[{p2,p4,V2,V3,T2,T3,T4,pV1,pV2,pV3,pV4,adiKomp,adiExp,plot,p,V},
V2=V1/ϵ;V3=V2*ϕ;
p2=p1*(ϵ^κ);p4=p2*(V3/V1)^κ;

adiKomp=p1*(V1^κ);adiExp=p2*(V3^κ);
pV1={V1,p1};
pV2={V2,p2};
pV3={V3,p2};
pV4={V1,p4};
plot=Table[{V,adiKomp*V^(-κ)},{V,V2,V1,.001}];

Which[
time≤10,V=V1-time*(((V1-V2)/10));p=adiKomp*V^(-κ),
time≤20,p=p2;V=V2+(time-10)*(((V3-V2)/10)),
time≤30,V=V3+(time-20)*(((V1-V3)/10));p=adiExp*V^(-κ),
time≤40,p=p4-(time-30)*(((p4-p1)/10));V=V1,
_,None];

Show[{
Plot[adiKomp*V^(-κ),{V,V1,V2},PlotStyle->{Blue,Thick}],
Plot[adiExp*V^(-κ),{V,V3,V1},PlotStyle->{Blue,Thick}],
Graphics[{

Blue,Thick,Line[{{pV2,pV3},{pV4,pV1}}],

Red,PointSize[.05],Point[{V,p}]
}]
},AxesOrigin->{0,0},PlotRange->All,AspectRatio->1,PlotLabel->"p-V diagram",ImageSize->imgS]
]

dieselPT[p1_,V1_,T1_,ϵ_,ϕ_,κ_]:=
Module[{p2,p4,V2,V3,T2,T3,T4,pV1,pV2,pV3,pV4,adiKomp,adiExp,isoch},
V2=V1/ϵ;V3=V2*ϕ;
p2=p1*(ϵ^κ);
p4=p2*(V3/V1)^κ;
T2=T1*ϵ^(κ-1);
T3=T2*ϕ;
T4=T3*(V3/V1)^(κ-1);
adiKomp=(p1^(1-κ))*(T1^κ);adiExp=(p2^(1-κ))*(T3^κ);isoch=(p4*V1)/T4;
pV1={T1,p1};
pV2={T2,p2};
pV3={T3,p2};
pV4={T4,p4};

Show[{
Plot[((adiKomp/(T^κ))^(1/((1-κ)))),{T,T1,T2},PlotStyle->{Blue,Thick}],
Plot[((adiExp/(T^κ))^(1/((1-κ)))),{T,T3,T4},PlotStyle->{Blue,Thick}],
Plot[(T*isoch)/V1,{T,T1,T4},PlotStyle->{Blue,Thick}],
Graphics[{

Blue,Thick,Line[{pV2,pV3}]
}]
},AxesOrigin->{0,0},PlotRange->All,AspectRatio->1,Axes->True,PlotLabel->"p-T diagram"]
]

dieselPTTime[p1_,V1_,T1_,ϵ_,ϕ_,κ_,time_,imgS_]:=
Module[{p2,p4,V2,V3,T2,T3,T4,pV1,pV2,pV3,pV4,adiKomp,adiExp,isoch,p,V,T},
V2=V1/ϵ;V3=V2*ϕ;
p2=p1*(ϵ^κ);
p4=p2*(V3/V1)^κ;
T2=T1*ϵ^(κ-1);
T3=T2*ϕ;
T4=T3*(V3/V1)^(κ-1);
adiKomp=(p1^(1-κ))*(T1^κ);adiExp=(p2^(1-κ))*(T3^κ);isoch=(p4*V1)/T4;
pV1={T1,p1};
pV2={T2,p2};
pV3={T3,p2};
pV4={T4,p4};

Which[
time≤10,T=T1+time*(((T2-T1)/10));p=((adiKomp/(T^κ))^(1/((1-κ)))),
time≤20,p=p2;T=T2+(time-10)*(((T3-T2)/10)),
time≤30,T=T3-(time-20)*(((T3-T4)/10));p=((adiExp/(T^κ))^(1/((1-κ)))),
time≤40,p=p4-(time-30)*(((p4-p1)/10));V=V1;T=T4-(time-30)*(((T4-T1)/10)),
_,None];

Show[{
Plot[((adiKomp/(T^κ))^(1/((1-κ)))),{T,T1,T2},PlotStyle->{Blue,Thick}],
Plot[((adiExp/(T^κ))^(1/((1-κ)))),{T,T3,T4},PlotStyle->{Blue,Thick}],
Plot[(T*isoch)/V1,{T,T1,T4},PlotStyle->{Blue,Thick}],
Graphics[{

Blue,Thick,Line[{pV2,pV3}],

Red,PointSize[.05],Point[{T,p}]
}]
},AxesOrigin->{0,0},PlotRange->All,AspectRatio->1,Axes->True,PlotLabel->"p-T diagram",ImageSize->imgS]
]

dieselTS[p1_,V1_,T1_,s1_,ϵ_,ϕ_,κ_,cp_]:=
Module[{p2,p4,V2,V3,T2,T3,T4,S1=s1,S2,S3,S4,TS1,TS2,TS3,TS4,adiKomp,adiExp,isoch,isob},
V2=V1/ϵ;V3=V2*ϕ;
p2=p1*(ϵ^κ);
p4=p2*(V3/V1)^κ;
T2=T1*ϵ^(κ-1);
T3=T2*ϕ;
T4=T3*(V3/V1)^(κ-1);
isob=Table[{cp*Log[tempr/T2],tempr},{tempr,T2,T3}];
isoch=Table[{(cp/κ)*Log[tempr/T1],tempr},{tempr,T1,T4}];
adiKomp=(p1^(1-κ))*(T1^κ);adiExp=(p2^(1-κ))*(T3^κ);

TS1={S1,T1};
TS2={S1,T2};
TS3={S1+cp*Log[T3/T2],T3};
TS4={S1+(cp/κ)*Log[T4/T1],T4};

Show[{
ListLinePlot[isoch,PlotStyle->{Blue,Thick}],
ListLinePlot[isob,PlotStyle->{Blue,Thick}],
Graphics[{

Blue,Thick,Line[{{TS1,TS2},{TS3,TS4}}]
}]
},AxesOrigin->{-1,0},PlotRange->All,AspectRatio->1,Axes->True,PlotLabel->"T-s diagram"]
]

dieselTSTime[p1_,V1_,T1_,s1_,ϵ_,ϕ_,κ_,cp_,time_,imgS_]:=
Module[{p2,p4,V2,V3,T2,T3,T4,S1=s1,S2,S3,S4,TS1,TS2,TS3,TS4,adiKomp,adiExp,isoch,isob,T,s},
V2=V1/ϵ;V3=V2*ϕ;
p2=p1*(ϵ^κ);
p4=p2*(V3/V1)^κ;
T2=T1*ϵ^(κ-1);
T3=T2*ϕ;
T4=T3*(V3/V1)^(κ-1);
S3=cp*Log[T3/T2];S4=(cp/κ)*Log[T4/T1];
isob=Table[{cp*Log[tempr/T2],tempr},{tempr,T2,T3}];
isoch=Table[{(cp/κ)*Log[tempr/T1],tempr},{tempr,T1,T4}];
adiKomp=(p1^(1-κ))*(T1^κ);adiExp=(p2^(1-κ))*(T3^κ);

TS1={S1,T1};
TS2={S1,T2};
TS3={S3,T3};
TS4={S4,T4};

Which[
time≤10,s=S1;T=T1+time*(((T2-T1)/10)),
time≤20,T=T2+(time-10)*(((T3-T2)/10));s=S1+cp*Log[T/T2],
time≤30,s=S3;T=T3-(time-20)*(((T3-T4)/10)),
time≤40,T=T4-(time-30)*(((T4-T1)/10));s=S1+(cp/κ)*Log[T/T1],
_,None];

Show[{
ListLinePlot[isoch,PlotStyle->{Blue,Thick}],
ListLinePlot[isob,PlotStyle->{Blue,Thick}],
Graphics[{

Blue,Thick,Line[{{TS1,TS2},{TS3,TS4}}],

Red,PointSize[.05],Point[{s,T}]
}]
},AxesOrigin->{-1,0},PlotRange->All,AspectRatio->1,Axes->True,PlotLabel->"T-s diagram",ImageSize->imgS]
]

pistCiferniky[p1_,V1_,T1_,ϵ_,ϕ_,κ_,time_,imgS_]:=
Module[{p2,p4,V2,V3,T2,T3,T4,pV1,pV2,pV3,pV4,adiKomp,adiExp,p,V,T},
V2=V1/ϵ;V3=V2*ϕ;
p2=p1*(ϵ^κ);p4=p2*(V3/V1)^κ;
T2=T1*ϵ^(κ-1);
T3=T2*ϕ;
T4=T3*(V3/V1)^(κ-1);

pV1={V1,p1};
pV2={V2,p2};
pV3={V3,p2};
pV4={V1,p4};
Which[
time≤10,p=p1+time*(((p2-p1)/10));V=V1-time*(((V1-V2)/10));T=T1+time*(((T2-T1)/10)),
time≤20,p=p2;V=V2+(time-10)*(((V3-V2)/10));T=T2+(time-10)*(((T3-T2)/10)),
time≤30,p=p2-(time-20)*(((p2-p4)/10));V=V3+(time-20)*(((V1-V3)/10));T=T3-(time-20)*(((T3-T4)/10)),
time≤40,p=p4-(time-30)*(((p4-p1)/10));V=V1;T=T4-(time-30)*(((T4-T1)/10)),
_,None];
GraphicsRow[{AngularGauge[V,{V2,V1},ImageSize->250,GaugeLabels->Style["V [(cm^3)]",Bold,14]],AngularGauge[p,{p1,p2},ImageSize->250,GaugeLabels -> Style["p [kPa]",Bold,14]],ThermometerGauge[T,{T1,T3},ImageSize->Medium,GaugeLabels ->Style[ "T [K]",Bold,14]]},ImageSize->imgS]
]

valec[p1_,V1_,ϵ_,ϕ_,κ_,time_]:=
Module[
{p2,p4,V2,V3,T2,T3,T4,pV1,pV2,pV3,pV4,adiKomp,adiExp,p,V,
valecVyska=8.5,tlIn=3,tlOut=.5,valecObjem=6,moveObjem,pistVyska=2,hridelDelka=2,hridelS,hridelR},

V2=V1/ϵ;V3=V2*ϕ;
p2=p1*(ϵ^κ);p4=p2*(V3/V1)^κ;
adiKomp=p1*(V1^κ);adiExp=p2*(V3^κ);
pV1={V1,p1};
pV2={V2,p2};
pV3={V3,p2};
pV4={V1,p4};
Which[
time≤10,p=p1+time*(((p2-p1)/10));V=V1-time*(((V1-V2)/10)),
time≤20,p=p2;V=V2+(time-10)*(((V3-V2)/10)),
time≤30,p=p2-(time-20)*(((p2-p4)/10));V=V3+(time-20)*(((V1-V3)/10)),
time≤40,p=p4-(time-30)*(((p4-p1)/10));V=V1,
_,None];
moveObjem=valecObjem*(V/V1);

Graphics3D[{
{FaceForm[Opacity[0.3]],Cylinder[{{0,0,0},{0,0,valecVyska}},tlIn+tlOut]},
{FaceForm[Opacity[0.3]],Cylinder[{{0,0,tlOut},{0,0,(valecVyska)}},tlIn]},
Cylinder[{{0,0,(moveObjem+tlOut)},{0,0,(tlOut+moveObjem+pistVyska)}},tlIn]

}]
]

valecAhridel[p1_,V1_,ϵ_,ϕ_,κ_,time_]:=
Module[
{p2,p4,V2,V3,T2,T3,T4,pV1,pV2,pV3,pV4,adiKomp,adiExp,p,V,
valecVyska=8.5,tlIn=3,tlOut=.5,valecObjem=6,moveObjem,pistVyska=2,
hridelDelka,hridelS,hridelR},

V2=V1/ϵ;V3=V2*ϕ;
p2=p1*(ϵ^κ);p4=p2*(V3/V1)^κ;
adiKomp=p1*(V1^κ);adiExp=p2*(V3^κ);
pV1={V1,p1};
pV2={V2,p2};
pV3={V3,p2};
pV4={V1,p4};
Which[
time≤18,p=p1+time*(((p2-p1)/18));V=V1-time*(((V1-V2)/18)),
time≤20,p=p2;V=V2+(time-18)*(((V3-V2)/2)),
time≤36,p=p2-(time-20)*(((p2-p4)/16));V=V3+(time-20)*(((V1-V3)/16)),
time≤40,p=p4-(time-36)*(((p4-p1)/4));V=V1,
_,None];
moveObjem=valecObjem*(V/V1);

hridelR=((valecObjem(ϵ-1))/(2ϵ));
hridelDelka=6;
hridelS=tlOut+valecObjem+pistVyska+hridelDelka-hridelR;

Graphics3D[{
{FaceForm[Opacity[0.3]],Cylinder[{{0,0,0},{0,0,valecVyska}},tlIn+tlOut]},
{FaceForm[Opacity[0.15]],Cylinder[{{0,0,tlOut},{0,0,(valecVyska)}},tlIn]},
Cylinder[{{0,0,(moveObjem+tlOut)},{0,0,(tlOut+moveObjem+pistVyska)}},tlIn],

Cylinder[{{-1,0,hridelS},{-.1,0,hridelS}},hridelR+.2],
PointSize[Large],Red,Point[{0,0,tlOut+pistVyska+moveObjem}],
{Blue,Point[{0,Bx,By}],
Gray,Tube[{{0,0,tlOut+pistVyska+moveObjem},{0,Bx,By}},.1],
Cylinder[{{-.1,Bx,By},{.1,Bx,By}},.15]}/.If[time≤10,
First@Solve[{Sqrt[(0-Bx)^2+(moveObjem+tlOut+pistVyska-By)^2]==hridelDelka,Sqrt[(0-Bx)^2+(hridelS-By)^2]==hridelR},{Bx,By}],
Last@Solve[{Sqrt[(0-Bx)^2+(moveObjem+tlOut+pistVyska-By)^2]==hridelDelka,Sqrt[(0-Bx)^2+(hridelS-By)^2]==hridelR},{Bx,By}]
]
},ImageSize->Large,Axes->True]
]